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We use 3-D radiative transfer models to show the effects of clumpy circum- 
stellar material on the observed infrared colors of high mass stars embedded in 
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molecular clouds. We highlight differences between 3-D clumpy and 1-D smooth 
models which can affect the interpretation of data. We discuss several important 
properties of the emergent spectral energy distribution (SED): More near-infrared 
light (scattered and direct from the central source) can escape than in smooth 
^ ■ 1-D models. The near- and mid-infrared SED of the same object can vary signifi- 

cantly with viewing angle, depending on the clump geometry along the sightline. 
Q\ ■ Even the wavelength-integrated flux can vary with angle by more than a factor 
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of two. Objects with the same average circumstellar dust distribution can have 
very different near- and mid-IR SEDs depending on the clump geometry and the 
proximity of the most massive clump to the central source. 

Although dumpiness can cause similar objects to have very different SEDs, 
there are some observable trends. Near- and mid-infrared colors are sensi- 
tive to the weighted average distance of clumps from the central source and 
to the magnitude of clumpy density variations (smooth-to-clumpy ratio). Far- 
infrared emission remains a robust measure of the total dust mass. We present 
simulated SEDs, colors, and images for 2MASS and Spitzer filters. We com- 
pare to observations of some UCH II regions and find that 3-D clumpy mod- 
els fit better than smooth models. In particular, clumpy models with fractal 
dimensions in the range 2.3-2.8, smooth to clumpy ratios of <50%, and den- 
sity distributions with shallow average radial density profiles fit the SEDs best 
{(p) ocr Q ,-1.0 < a < 0.0). 
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1. Introduction 

The important question of how stars form is frequently investigated with infrared imag- 
ing, for both astrophysical and practical reasons. Young stellar objects (YSOs) are usually 
still deeply embedded in their natal molecular material. The distribution of this gas and dust 
(e.g. in an accretion disk) can reveal a great deal about the formation process, and this hot 
dust close to the star emits primarily in the near and mid-infrared (for this paper, we define 
near-IR = 1< A <3/im and mid-IR =3< A <30/im). Although the bulk of the dust mass is 
often at cooler temperatures, with SEDs peaking in the far-infrared (~100/im), probing the 
warmer dust is both advantageous and potentially dangerous because the NIR/MIR SED 
is more sensitive to the dust distribution and geometric effects. The infrared is observa- 
tionally advantageous because at shorter wavelengths, YSOs are completely invisible, and at 
longer wavelengths, detector technology is less mature and the spatial resolution of (typically 
diffraction-limited) observatories begins to suffer. Many current and recent observatories are 
infrared-optimized, such as ISO and ground-based telescopes with adaptive optics. Particu- 
larly exciting is the successful operation of the Spitzer Space Telescope, which is producing a 
wealth of MIR data on star formation due to its unprecedented mapping efficiency (Allen et 
al. 2004; Megeath et al. 2004; Reach et al. 2004; Whitney et al. 2004a). The IRAC instru- 
ment (Fazio et al. 2004) provides simultaneous < 2" imaging at 3.6, 4.5, 5.8, and 8.0 /xm, 
and MIPS (Rieke et al. 2004) provides 6"images at 24/im, and longer wavelengths at lower 
resolution. It is critical that we understand how circumstellar dust geometry can affect our 
interpretation of all of this NIR/MIR data. 

In a previous set of papers (Whitney et al. 2003a,b, 2004b) we showed how 2-D geome- 
tries (disks, rotationally flattened envelopes, bipolar cavities) affect the near-IR and mid-IR 
colors of Young Stellar Objects (YSOs) of various stellar masses and evolutionary states. 
Those models show that there is not always a direct correspondence of the 1-10 /im colors 
or spectral index with evolutionary state, due to several effects such as scattered light, incli- 
nation, aperture size, and stellar temperature. There is a trend for younger sources to have 
redder colors, but there is substantial overlap and extreme exceptions to the trends. Sepa- 
ration of these populations in a color-color diagram generally requires supplemental longer 
(A > 20/im) wavelength observations. If only 1 — 10/im SEDs are available, for example, in 
the Galactic Legacy Infrared Midplane Extraordinaire (GLIMPSE, a Spitzer project; Ben- 
jamin et al. 2003), model fits to the spectral energy distributions (SEDs) provide a more 
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promising avenue than color-color diagrams for interpreting the SEDs since they use the 
most information simultaneously (Robitaille et al. 2005, in prep.). 

This paper focuses on UltraCompact H II regions, or the later stages of high mass 
(O and early B) star formation; although we have not yet explored the issue, we expect 
that the dumpiness of envelopes is important also in younger high-mass objects, and per- 
haps in somewhat lower-mass protostars as well. It is believed that the central sources of 
UltraCompactH II regions (UCH lis; Churchwell et. al 2004) have probably reached the 
main sequence, possibly halted accretion and dispersed infalling material, but remain deeply 
embedded in their natal molecular cloud. These stars are so luminous that they heat up 
large volumes (~ 3 pc radius) and masses of dust and gas, which reprocess the stellar ra- 
diation towards longer wavelengths. Observations of high mass YSOs and UCH II regions 
are usually modeled using 1-D radiative transfer codes (Wolfire & Churchwell 1994; Faison 
et al. 1998; Hatchell et al. 2000; van der Tak et al. 2000; Beuther et al. 2002; Mueller et al. 
2002). These are adequate for modeling long wavelength observations (> lOOyum), but at 
shorter wavelengths, the inhomogeneous nature of the surrounding interstellar medium (e.g. 
Elmegreen & Falgarone 1996) likely has a large effect on the SED. This is demonstrated 
in the attempts by Faison et al. (1998) (hereafter F98) to fit the mid-IR spectroscopy of 
several UCH II sources with 1-D smooth models. For a given envelope mass set by the long 
wavelength flux, the models consistently produced too deep a 10 fxm silicate feature and 
too little short wavelength flux, even after decreasing the silicate abundances by a factor of 
two. Changing the radial density profile is insufficient to bring the models into agreement 
with the data. F98 use constant-density spheres with large inner holes and obtain a poor 
match. Some of the other models cited above use negative power-law density gradients to 
fit spatially resolved intensity profiles; such 1-D models still have difficulties fitting UCH II 
SEDs over a wide wavelength range (§3.2). 

In §2 we describe our models. Three-dimensional nonhomogeneous models differ signif- 
icantly from 1-D smooth models with the same amount of circumstellar dust. Some authors 
have computed SEDs in two-phase clumpy distributions (Misselt et al. 2001; Wolf et al. 
1998). In an attempt to calculate as realistic models as possible, in this paper we calcu- 
late SEDs in hierarchically clumped density structures, following a procedure for calculating 
fractal density described by Elmegreen (1997) and Mathis, Whitney, & Wood (2002) and 
using a volume fractal dimension D similar to that calculated for the interstellar medium 
(D = 2.3 — 2.6; Elmegreen & Falgarone 1996). In §3.3 we describe the properties of the 3-D 
models in terms of synthetic observations: high-resolution images, spectral energy distribu- 
tions, and photometric fluxes analyzed using color-color and color-magnitude diagrams. We 
discuss how the properties vary with viewing angle (§3.4), and show in detail how the SEDs 
of 1-D models cannot reproduce the characteristics of 3-D models (§3.6). In §4.1 we show 
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that despite the variations, there are some trends that allow observations to probe physical 
conditions. In §5 we demonstrate our findings on real data, and give conclusions in §6. 

2. Model Construction 

We use a Monte Carlo Radiative transfer code similar to that described in Whitney et al. 
(2003a) (Paper I), but using a 3-D grid and incorporating a hierarchically clumped density 
distribution. We use the radiative equilibrium method developed by Bjorkman & Wood 
(2001). The models solve for the 3-D temperature distribution, conserve flux absolutely, and 
accurately compute scattering and polarization using arbitrary scattering phase functions. 
The output radiation field is binned into 200 angular directions (10 polar intervals and 20 
azimuth intervals), i.e. each model produces 200 SEDs corresponding to different viewing 
angles. We do not currently include emission from polycyclic aromatic hydrocarbons (PAHs) 
or transiently heated grains; these likely affect fluxes in the IRAC wavelength range (3-8/zm) . 
The only source of calculation error in our models is photon counting statistics - running 
more photons produces higher signal-to-noise spectra. The models produced for figures in 
this paper took 12-24 hours each to run (on a 2 GHz PC and code compiled with g77). Draft 
quality results sufficient to explore trends take <2 hours. 

We specify circumstellar dust distributions with a clumpy and smooth component, with 
the relative contribution ranging from zero (very clumpy with evacuated holes) to one (con- 
stant density). The clump distribution is specified using the algorithm of Elmegreen (1997) 
as implemented in Mathis, Whitney, & Wood (2002): N mass particles are randomly placed 
in a simulation grid of size L. N more particles are randomly placed, each within a distance 
L/(2A) of one of the particles in the previous round (A is a parameter which determines 
the approximate fractal dimension of the resulting distribution). This process is repeated 
twice more, so that iV 4 particles are placed. The mass distribution is not a true fractal (self- 
similar on all scales), but merely hierarchically clumped (self-similar over about an order of 
magnitude in scale). Nevertheless, the quantity D = log(iV)/ log(A) is essentially the fractal 
dimension. The precise value of D does not dramatically affect any of our results - this will 
be demonstrated in §4.1.1, after discussing more universal properties of the models. 

Different grain models could be used in different parts of the dust cocoon (e.g. larger 
grains in dense clumps), as in the young stellar objects modeled in Whitney et al. (2003b). 
For simplicity, we use the same grain model throughout the envelope; as described in Pa- 
per I, the grain distribution was fit to the extinction curve typical of molecular clouds with 
-Ry=4.3 using a maximum entropy method similar to Kim, Martin, & Hendry (1994). We 
use spherical grains with dielectric functions for silicate and graphite from Laor & Draine 
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(1993). The dust sublimation temperature is chosen to be 1600 K. 

For simplicity, all of the clumpy models shown in this paper are for a massive star, with 
the same amount of circumstellar dust, and the same photospheric spectrum of the central 
source (the appropriate Kurucz model). The results presented are general to all massive 
young stars, but we wish to isolate the effects of differing clumpy circumstellar distribution 
from the effects of different stellar masses and temperatures. The latter effect is discussed in 
(Whitney et al. 2003b) for 2-D circumstellar distributions, and further discussion of clumpy 
dust in the context of embedded protoclusters is presented in Johnson et al 2005 (submitted). 

In order to provide a basis for comparison, we ran numerous models with the smooth-to- 
clumpy density ratio equal to one, i.e. const ant- density spheres. In particular, we compare 
our results with the popular DUSTY code (Ivezic & Elitzur 1997) in Figure 1. Across a 
range of parameters, the SEDs are generally indistinguishable between the two codes, and 
quantitatively agree to a few percent. Our code properly calculates nonisotropic scattering 
(Paper 1), whereas DUSTY uses an isotropic scattering function. For models with low 
(ry < 30) optical depth, this can increase the differences between the two codes to 10-15%, 
but if we run our code with isotropic scattering the results are then indistinguishable. We 
will demonstrate that differences between 1-D and 3-D models are far greater than any of 
these subtle (few percent) effects. 

3. Results of a Canonical Clumpy Model 

3.1. Density and Temperature Distributions 

Figure 2 shows an azimuthal slice through the density grid of our canonical model, a 
main sequence 06-07 star (T*=41000 K, L 7t =2.5xlO 5 L ) in a dust cocoon extending from 
0.0001 pc to 2.5 pc, with a smooth-to-clumpy ratio of 0.1 and fractal dimension 2.6. The 
cocoon is representative of the inferred properties of massive YSOs: a thousand solar masses 
of dust in about a 1/3 parsec region. (At typical Galactic distances of 5kpc, l/3pc=13"and 
2.5pc=100".) As massive YSOs can heat a large part of their natal clouds, we chose to 
include the larger envelope out to 2.5pc and 50000 M Q ; this does not qualitatively affect 
our conclusions. The visual extinction averaged over all sightlines is A V =IS1, but varies 
between 13 to 401 magnitudes depending on viewing angle. The voids at small radii allow 
stellar radiation to propagate further, heating the inner faces of clumps. This can result 
in a greater mass of warm dust than if all of the stellar optical and UV radiation was 
intercepted near the star; in other words, the Ty — 1 surface is irregular and will enclose a 
different amount of mass than in a spherically symmetric model. This is clearly visible in 
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Fig. 1. — The first panel shows a 1-D model (spherically symmetric) calculated with our code 
(solid line) and DUSTY (dashed line). The source is a main sequence 06-07 star (T*=41000 K) 
in a 50000 M Q dust cocoon extending from 0.0001 to 2.5 pc. The optical depth of the cocoon is 
Ty=117.5, and the temperature on the inner radius is 700 K. The second panel shows the fractional 
difference between our model and one calculated with DUSTY. Variations of a few percent are 
mostly due to noise in our spectrum. 

the temperature map, which shows cool dust in the shadows of clumps, and heated inner 
surfaces. 

The radial temperature profile is sometimes used as a diagnostic when young stellar ob- 
jects are modeled. As one would expect, the radial temperature profiles for clumpy models 
vary greatly depending on the placement of the clumps of dust. Even the radial temper- 
ature profiles averaged over all angles show significant variation with the particular clump 
geometry, and differ from the profile of a smooth (1-D) model with the same amount of 
circumstellar mass. Figure 3 illustrates this comparison with smooth and clumpy models 
that have the same amount of circumstellar mass. The smooth model's temperature gradient 
is shallower than the clumpy models in the inner regions, because the clumpy models tend 
to have lower average density there. However, a smooth model with a completely evacu- 
ated hole out to the distance of the nearest dense clump in the clumpy models also shows 
a faster temperature dropoff than the clumpy models. This is the effect mentioned above, 
that the r=l surface for intercepting stellar radiation is irregular, and clump surfaces can 
be heated at larger distances from the star than in a spherically symmetric model with the 
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Fig. 2. — Density and temperature distributions in a slice of a typical model. Both quantities 
are shown on a logarithmic scale, with the ranges shown by the color bars. The spatial scale is in 
parsecs. Two particular sightlines discussed in the text ("A" and "B", §3.4) are shown as dotted 
lines (although we will discuss a sightline as if it were a single viewing angle, radiation leaving the 
simulation is binned over the finite range of angles indicated by the solid lines). 
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most similar radial dust distribution. For the dust used here, temperature is expected to 
drop off as r -0 33 or r~ 0A in the optically thin and thick smooth cases, respectively (more 
discussion in Whitney et al. 2003a). The different slopes for different clumpy models likely 
reflect transitions between optically thin interclump material and optically thick clumps. 

3.2. SEDs 

We show here a comparison of SEDs of the clumpy and smooth model for our canonical 
parameter set discussed in § 2. As Figure 4 shows, our clumpy model differs from the smooth 
model most significantly in the NIR and MIR. The voids between the clumps (Figure 2) allow 
NIR flux to scatter through and escape the envelope. MIR flux also can escape along certain 
sightlines. The 10 /zm silicate feature varies greatly between the different sightlines of the 
clumpy model, but the average depth is lower than that of the smooth model. The slope of 
the SED in the NIR and MIR is flatter on average in the clumpy model than the smooth. 
As we will show in later sections, these two features of the clumpy models agree better with 
observational data of several UCH II sources. 

Figure 2 showed a slice of the density and temperature structure for a typical clumpy 
model, with two particular sightlines marked. The SEDs for all 200 sightlines in this model 
are shown in Figure 4 (in gray), with the two sightlines highlighted. The output spectra 
are binned over a finite range in angle, indicated in Fig. 2 for each sightline by two solid 
lines bracketing the dotted line. It is important to remember that although we have only 
shown the sightline intersecting the central source, the SED at each viewing angle is due 
to emission in that direction from the entire object, i.e. all rays parallel to the range of 
angles indicated. The SED with lesser near-infrared flux corresponds to the lower ("A") of 
the two sightlines in Fig. 2 - there are clumps in the outer envelope in this direction that 
obscure the inner regions (this cooler dust also causes the lO/im silicate absorption present 
in the SED, Fig. 4). Along the upper sightline "B", there are few clouds at a radius greater 
than 0.4 pc, and much more NIR radiation can scatter out from the inner regions (dashed 
SED, Fig. 4). The orientation of the two sightlines is particular to this model - there is no 
preferred direction to the clump distribution in general. 

3.3. Synthetic Images 

Figure 5 shows 3-color images of our canonical clumpy model as it would be imaged with 
a standard ground-based NIR camera, and with the IRAC and MIPS on Spitzer. Specif- 
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Fig. 3. — Comparison of spherically-averaged radial temperature and density profiles for smooth 
and clumpy models. All models shown here have the same amount of circumstellar dust within 
the same 2.5 pc radius. Two smooth models (solid and dashed lines) and 4 clumpy models (long 
dashed, dot-dashed, dotted, and triple-dot-dashed lines) are shown; the density and temperature 
profiles for the clumpy models are averages over all angles. The upper plot shows the temperature 
profiles, and the lower the corresponding density profiles. The smooth model which extends all the 
way to the dust destruction radius at constant density (dashed line) has the shallowest gradient, 
because it has the highest density in the inner regions of the four clumpy models. The temperature 
of a smooth model with a completely evacuated inner hole (solid line) also falls off faster with radius 
than any clumpy model, because the inner faces of clumps at varying distances are illuminated and 
heated, whereas all stellar radiation of a given wavelength is intercepted at the same radius in a 
smooth model. 
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wavelength [/tm] 

Fig. 4. — SEDs of clumpy and smooth models. The first panel shows the range of SEDs for 
different sightlines of the canonical clumpy model. The grey scale is darker where more similar 
SEDs coincide. SEDs of two particular sightlines ("A" and "B", see Figs. 2 and 5) are highlighted. 
The second panel compares the clumpy model to a smooth model (dashed line) with the same dust 
mass and inner and outer radii. The average sightline of the clumpy model is shown as a solid 
line. The third panel makes the same comparison between a clumpy and smooth model, but now 
with an average radial density power-law gradient (p oc r~ l ). For simplicity we mostly discuss the 
differences between smooth and clumpy models with no radial dependence of the average density. 
In the model with a gradient, the differences are even larger - both models have higher average 
optical depth (than the models without a radial gradient) and more hot dust close to the star 
emitting in the NIR, but that emission only escapes and is observable in the clumpy model, which 
has low-density sightlines. In the smooth model all of that extra emission is preprocessed further 
out in the envelope. 
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ically, the images are constructed by convolving the simulated output radiation field with 
the (J,H,K S ), ([3.6], [4.5], [8.0]), and ([24], [70], [160]) filters. The top row shows a spherically 
symmetric model, and the second two rows are a clumpy model viewed from different angles. 
Clearly, the greatest difference is in the near-IR - the clumpy model is >6 orders of magni- 
tude brighter at these wavelengths. The clumpy model is also brighter in the mid-infrared 
(IRAC wavelengths), but more striking is the very different color - the smooth model is 
much redder, with higher average extinction than the clumpy model, and cooler dust farther 
out in the cocoon on average. Differences are less significant at longer wavelengths, which 
mostly trace the total dust mass, and are less sensitive to its distribution. This agrees with 
the SEDs shown in Figure 4. 



3.4. Integrated SED 

Not only does the shape of the spectral energy distribution change with viewing angle, 
but the total wavelength-integrated flux can change significantly as well. Figure 6 shows 
contours of the bolometric flux as a function of <fi and 9, the azimuthal and polar/altitude 
viewing angles, for a typical clumpy model. The luminosity that an observer would calculate 
by (incorrectly) assuming isotropic emission varies from 60% of the true luminosity along 
some sightlines to nearly twice the true value along other sightlines. The "dim" sightlines 
are those along which a dense clump is located at a moderate to large distance from the 
central source. Dense clumps along one's sightline, but very close to the central source, 
do not necessarily result in an overall dimming, because short-wavelength light can scatter 
around the clump, the viewer may still be able to see uneclipsed hot dust in the cocoon. 
We note that this variation of the bolometric flux was also found in the low-mass protostar 
model of Whitney et al. (2003a). 



3.5. Colors 

Observationally, the SED of young stellar objects is usually only sampled with a few 
broad-band filters, so it is important to show how the variation apparent in Figure 4 trans- 
lates into a range of photometric colors. With the successful launch of Spitzer, the availability 
of 2MASS, and widespread access to deep near-infrared imaging, photometry of a great num- 
ber of young stellar objects is available in the 10 filters J, H, K, IRAC [3.6], [4.5], [5.8], and 
MIPS [24], [70], [160], so we concentrate on these (though we note that in many cases the 
160 fim band is not very observationally useful because of its low resolution and saturation 
limit). For conversion of fluxes to magnitudes we use 1594., 1083., and 666.7 Jy for 2MASS 



Fig. 5. — Synthetic 3-color images of a massive star embedded in clumpy cloud. The top row is 
a spherically symmetric model, and the second two rows are two views of a clumpy model with 
the same mass of circumstellar dust. (The two sightlines are those marked in Figure 2.) The NIR 
(left panels) images are composed of the emergent spectrum convolved with the 2MASS J,H, and 
K s filters as blue, red, and green respectively and a logarithmic stretch. The MIR (center panels) 
images are similarly composed from synthetic IRAC [3.6], [4.5], and [8.0] images, and the FIR 
(right panels) images from synthetic MIPS [24], [70], and [160] images. We have not convolved 
these images with Spitzer resolution (1.5"at 6/iin and diffraction limited at higher wavelengths); 
the effect would be negligible for the NIR and MIR image, but would begin to strongly affect the 
FIR image at the assumed distance of 5 kpc. 
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Fig. 6. — Observed bolometric flux as a function of polar /altitude (9) and azimuthal ((f)) viewing 
angle. The light spots are angles along which the observed luminosity (if the observer assumes 
isotropic emission) is twice the true value, and in the dark regions the observed value is 60% of the 
true value. 
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J, H, and K s , 277.5, 179.5, 116.6, 63.1, and 7.16 Jy for the four IRAC bands and MIPS [24] 
(M. Cohen 2004, private communication), and an estimated 0.75 and 0.2 Jy for MIPS [70] 
and [160]. We also calculate the IRAS fluxes of our models and discuss these briefly in §5. 

Figure 7 shows two color-color diagrams for our example typical clumpy model, viewed 
from 200 different angles. The near-infrared bands vary the most with viewing angle, but 
colors formed with an IRAC band and MIPS [24] are also quite sensitive to the circumstellar 
dust distribution. This is because the MIPS bands are at long enough wavelengths that they 
sample predominantly thermal emission of fairly cool dust; the flux density mostly depends 
on the total dust mass and is less affected by the particular distribution. We also show the 
color-color diagram formed exclusively from IRAC bands, since many young stellar objects 
in the inner Galaxy are being imaged with IRAC by the GLIMPSE survey (Benjamin et al. 
2003; Churchwell et. al 2004; Whitney et al. 2004a). Colors can vary with viewing angle by 
more than a magnitude for the same object. The corresponding smooth model (same dust 
mass and inner and outer radii) is shown on the two CCDs, indicating that even though 
the variation between sightlines in the clumpy model is large, the difference between any 
sightline and a smooth model is even larger. The smooth model is extremely embedded and 
nearly invisible in the NIR, so the [K]-[3.6] color is far off the right side of the first panel. 
Even the IRAC-only colors of the smooth model are different from the clumpy model by 
more than a magnitude. 



3.6. Fitting 3-D data with 1-D Models 

The spectral energy distributions of 3-D clumpy models are significantly different from 
those of 1-D smooth models. We attempt to fit the output of our 3-D models with 1- 
D models, revealing several specific differences: 1) The derived physical parameters (dust 
mass, etc) are often very different, 2) the derived physical parameters for the same object 
depend on viewing angle, and 3) sometimes it is impossible to fit the (synthetic) data with 
a 1-D model at all. Figure 8 illustrates these points. The range of SEDs for a single clumpy 
model seen from different angles is shown as a gray scale, along with three 1-D models. The 
"true" physical parameters are 50000 M Q of circumstellar material around an O-type star 
(T* =41040 K, L* =2.54xl0 5 L ). Some very thin dust extends all the way in to 0.0001 pc 
(the density is never truly zero), and the closest clump (typical density of clumps is 10 5 cm -3 ) 
is at radius 0.035 pc, and temperature 250 K. The average optical depth to the central star is 
Ay =131, but that quantity ranges from 13 to 401 magnitudes depending on viewing angle. 

We have concentrated here on fitting the NIR and MIR parts of the SED, without 
heavily constraining the FIR. This is motivated by current observational capability, which 
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Fig. 7. — Color-color diagrams (CCDs) for different observations of a single modeled object, an 
0-type star with 50000 M of circumstellar material distributed in clumps between 0.0001 and 
2.5 pc radius (the same model as Figs. 2-6). The different points are for the 200 different viewing 
angles calculated for each model. In the first panel, colors across the broad wavelength range 
2-25//m are shown, and in the second panel colors derivable exclusively from IRAC photometry. 
The corresponding smooth model is shown as a filled triangle in the first panel, indicating that the 
actual point is far off the plot to the right, and as a filled square in the second panel. 
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Fig. 8. — Attempted 1-D fits to the SEDs of a 3-D clumpy model. The variation of SEDs with 
viewing angle is shown in grayscale as in Figure 4. The high and low extremes of the greyscale 
correspond to the two sightlines shown in Fig. 2 that have particularly high and low NIR flux 
densities ("B" and "A", respectively, see first panel of Fig 4). The dotted line is the best 1-D fit 
to SED "B". The dashed and dot-dashed are 1-D fits to SED "A". The 1-D model parameters do 
not reflect the true dust distribution very closely - see text for more discussion. 
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includes the widespread availability of NIR spectrographs and cameras and Spitzer's rela- 
tively high-resolution coverage of 3-30/im. IRAS and even the 70 and 160/im channels of 
MIPS on Spitzer are resolution-challenged, and it is more difficult to draw conclusions about 
massive protostars (generally at distances >1 kpc and in crowded environments) from these 
instruments. One can instead try to crudely fit the shape of the entire SED using a 1-D 
model; a reasonable match to the SED shape can be made for the least embedded sightlines 
as shown below. However, in general these 1-D fits either fail to reproduce the 10/im silicate 
feature in absorption or have MIR (~5/im) fluxes that are 1-3 orders of magnitude too low, 
and more seriously, the physical parameters derived are not the same as the 3-D model (size 
and mass of the cocoon). 

In order to quantify the differences between the SEDs produced by smooth 1-D and 
clumpy 3-D models, we constructed a grid of 5000 smooth models, varying the inner and 
outer radius, the optical depth or dust mass, and the radial density profile (power-law pro- 
files p{r) oc r a ). We then calculated the "x 2 " goodness-of-fit between each 1-D model and 
a selected 3-D sightline (and below, observation of a real YSO). We used an arbitrary "un- 
certainty" of 0.1 dex to calculate \ 2 f° r the model SEDs - the preferred parts of parameter 
space are not affected by this choice. Figure 9 shows this goodness of fit as a function of the 
optical depth ry, the inner and outer radii r in and r out , and the radial density power-law 
index a. 

The dotted line shows the best 1-D fit to the 3-D sightline with particularly high 
NIR/MIR flux density ("B" in Fig. 2). The best fitting models have fairly flat radial density 
profiles (a = 0.0; fit parameters are listed in Table 1). The preferred optical depths are 
much lower than the average optical depth of the 3-D model (21% for this best-fit model), 
so the dust mass would be severely underestimated. If one decreases the optical depth in the 
smooth model in order to increase the flux density at shorter wavelengths, one is also forced 
to increase the outer radius and total dust mass to maintain the observed FIR flux density. 
It is difficult to produce sensible 1-D models that have the 10/im silicate feature in emission 
as well as the shape of the SED at other wavelengths (in the 3-D models this is possible for 
viewing angles at which very hot dust near the star is not severely obscured.) 

It is even more difficult to fit the sightlines like "A" that have low MIR intensity using 
a 1-D model. The dot-dashed line shows the best fit constant-density model, which has a 
reasonably deep silicate absorption feature but not enough NIR flux and too much flux at 
20/im. This is a general feature with 1-D models that try to mimic typical sightlines of the 
3-D model: the optical depth must be kept moderate (similar or somewhat greater than the 
average value of the clumpy models) to produce the silicate absorption feature, and then 
there is no way to get emission shortward of 2/xm out of a reasonably embedded source 
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Fig. 9. — Goodness of fit of 1-D smooth models to two selected 3-D clumpy SEDs. The first two 
panels represent attempts to fit sightline "A" with a strong silicate absoprtion feature. The second 
two panels are for sightline "B" which shows silicate emission. In both cases, sightlines with fairly 
flat radial density profiles and larger inner and outer radii than the true model are preferred. The 
best fits have fairly poor "x 2 " values of 5 — 10. 
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Table 1. Parameters of 1-D models fitting 3-D models 



model 


density 






mass 


luminosity 


Average A v 


"x 2 " 




index a a 


[PC] 


[PC] 


M 


10 5 L 


to star 




clumpy 





0.0001 b 


2.5 


50000 


2.54 


131 




1-D fits to lower sightline 


"A" 














dashed line Fig 8 





0.047 


47. 


2.9xl0 7 


2.54 


217 


8.7 


dot- dash line Fig 8 


+1.5 


0.00039 


12. 


2.0xl0 6 


2.54 


163 


4.8 


luminosity unconstrained 





0.0005 


91 


5.6xl0 7 


0.6 


109 






-1 


0.0070 


7.0 


1.4xl0 5 


0.5 


200 




1-D fits to upper sightline 


"B" 














dotted line Fig 8 





0.0027 


2.7 


2.3xl0 5 


2.54 


27 


8.2 


luminosity unconstrained 





0.0023 


3.5 


31300 


3.5 


43 





a Smooth models were run with a power-law density gradient p oc r a . 

b Clumpy models have diffuse material extending in to the dust destruction radius, but the 
distance to the nearest dense (>10 5 cm~ 3 ) clump is farther out (0.035pc). 
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without an inhomogeneous envelope. In 1-D models the depth of the silicate feature is tied 
to the mass of the envelope, but three dimensions these are decoupled. 

The best-fitting models also have unreasonably large outer radii, and often, large inner 
holes or radial density profiles that increase with radius (a low-density inner region is similar 
to a hole of zero density). It is these last models (positive a) which produce the peak in the 
goodness-of-fit at very low r and very large r out - we do not consider these models physically 
reasonable, but plot one for completeness (dashed line). We note that 1-D fits to real data 
have the same characteristics - flat radial profiles and large outer and inner radii. (Our fits 
are discussed below, others such as F98 found the same result.) 

If one did not know the intrinsic luminosity of the source, that parameter would also 
likely be fit incorrectly. We did not construct a full grid of models with different luminosi- 
ties, but explored a range of parameter space fitting by eye, and found that the fit could 
be improved slightly over the SEDs shown in Figure 8 by changing the source luminosity 
(downward for sightline "A" and upward for sightline "B" - those model parameters are also 
listed in Table 1). We will see in later sections that 1-D models face similar difficulties in 
fitting real data that they do in fitting 3-D SEDs. 

4. More General Clumpy Models - Observable Trends 

4.1. Relating the SED and the dust distribution 

Although the same young stellar object can look different depending on the viewing 
angle, trends do exist that are strong enough to overcome that intrinsic scatter and probe 
the physical state of the circumstellar dust. Spatially averaged information about the density 
distribution and size of the inner cavity, the magnitude of density fluctuations relative to 
the mean, and the total dust mass are all observationally accessible. 

4-1.1. Trends with average Ay and inner cavity size 

One intuitive fact about circumstellar dust that remains true even in the presence of 
complex density distributions is that as the optical depth to the source of emission increases, 
the emergent near- and mid- infrared light decreases and the source appears redder. If one 
considers different sources with the same mass of circumstellar dust, then the size of an 
inner cavity or hole should be anticorrelated with the optical depth to the central source. 
In clumpy models, there is not a spherical inner hole, but the average distance of dust 
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from the central star varies depending on the specifics of the clump distribution. We can 
explore this effect in a natural way by running many models with different random seeds or 
realizations of the clump distribution. Figure 10 illustrates the results: We show 100 models 
with a smooth-to-clumpy ratio of 0.1 (mostly clumpy), and different random seeds of the 
clump distribution. The 20000 SEDs (200 viewing angles per model) vary widely at the 
shorter wavelengths, and the [5. 8] -[24] color shown in the lower panels can vary by over two 
magnitudes for a single model. The expected trend is nevertheless visible - models in which 
the clumps are largely further from the central source have somewhat lower average optical 
depth, and are somewhat bluer. The trend is quite clear in the angle-averaged properties of 
models (see the five highlighted models in particular), but the angle-averaged properties are 
not observable from Earth. 

The effective fractal dimension D offers another way to quantify the dust distribution - 
a larger effective dimension corresponds to more "fluffy" clouds with a larger filling factor. 
Most of our models use D=2.6, close to the effective dimension suggested by Elmegreen 
& Falgarone (1996) for the general interstellar medium, but we verify here that within a 
reasonable range, this choice does not have a strong effect on our results. Figure 11 shows 
slices of models with higher and lower D, as compared to the standard shown in Fig. 2. The 
SEDs and properties of models with a range of D are shown in Figure 12, where as before 
the full range of sightlines is shown and selected angle-averaged properties are highlighted. 
The variation of observable properties (SED, colors) is less than the variation as one changes 
the viewing angle for a single model. A fluffier or more filled cloud distribution results in 
somewhat redder colors and larger average optical depth, which is not surprising, but the 
effect is smaller than others examined in this paper, such as the cloud/intercloud density 
contrast (smooth-to-clumpy ratio). 
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Fig. 10. — One hundred clumpy models with the same amount of circumstellar mass in different 
random clump distributions. The upper panel shows the range of the 20000 sightlines (200 per 
model) as gray lines. Five particular models have been highlighted in color, and the angle-averaged 
SED plotted for those five. The lower panels show how the 5.8/xm-24/zm color varies as a function 
of the density weighted average radius (a measure of the distance from the central source at which 
most of the dust is located) and of the angle-averaged optical depth to the central source. The 200 
sightlines for each model are marked as dots in the lower panels, and the angle-averaged properties 
of the five highlighted models are shown as large diamonds. 




Fig. 11. — Models with effective fractal dimensions D=2A and 2.9 - compare to Fig. 2, which has 
D=2.6, the standard used in most of our models. Higher effective fractal dimension corresponds to 
"fluffier" clouds, a more filled cloud distribution, and lower density peaks. 
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Fig. 12. — The effective fractal dimension D does not have a strong effect on results (within the 
range shown, which is believed to be the relevant range for the ISM). The 200 sightlines for each 
of 7 models with different D are shown as gray scale in the upper panel, and dots in the lower 
panel. The angle-averaged properties are shown as thick SEDs and large diamonds. The variation 
in observable properties with D is less than the variation of a single model with viewing angle. 
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4.2. Trends with smooth-to-clumpy ratio 

To quantify the effects of a clumpy circumstellar medium in another way, we varied 
the ratio of the typical clump density to that of the smooth interclump medium. Figure 

13 shows the effect of varying the smooth-to-clumpy ratio for models the same particular 
random distribution of clouds. As in similar figures above, we show the range spanned by 
different viewing angles as well as the average properties to highlight the trend. As the 
smooth-to-clumpy ratio increases, the depth of the 10/im silicate feature increases, and the 
models become redder shortward of that feature. The MIR flux actually increases slightly for 
some parts of parameter space: the clumpy models can have a central hole with fairly thin 
material if there is not a clump immediately near the star, or a similar effect occurs if there 
is a clump at small radii but only one one side of the star. The effect of making such a model 
smoother is in fact to fill in that inner region with denser material, which can increase the 
mass of warm dust. In the smoother models, however, much of the radiation from such dust 
(at wavelengths shortward of a few microns) is extincted and reprocessed back in the far-IR 
part of the SED (see the rise in the peak at ~100/zm in Figure 13). The result is commonly 
the relative lowering of the 20/zm and NIR fluxes as shown in the Figure. Although there 
are large variations with viewing angle, the smoothness of the cocoon has a large enough 
effect that it can be observationally probed. To further illustrate that, we ran many models 
with different smooth-to-clumpy ratios and different random clump distributions. Figure 

14 shows the results. Particularly if one can observe a sample of similar objects, the mean 
color of the sample will reflect the smoothness of the circumstellar dust distribution - we 
will demonstrate in the next section. 



5. Comparison with data 

Now we apply the lessons of the previous section to real datasets, to determine what we 
can learn about the physical properties of ultracompact H II regions (embedded massive pro- 
tostars) using our clumpy models. Wood & Churchwell (1989) (hereafter WC89) examined 
the IRAS colors of ultracompact H II regions and determined that objects in a sufficiently 
red part of color space were likely to be massive protostars. This criterion has been used 
widely as the starting point for surveys of protostars and protostellar candidates, including 
thermal and maser molecular line targeted observations and imaging from radio continuum 
through the infrared. Although higher resolution more sensitive infrared data (MSX and 
Spitzer/ 'GLIMPSE) is beginning become available, IRAS is still a cornerstone dataset for 
studying dusty Galactic objects, so we calculate and describe the IRAS fluxes and colors of 
our models here. Figure 15 shows a color-color diagram similar to that of WC89. The WC89 
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Fig. 13. — SEDs for models with different ratios of smooth to clumpy density distribution, from 
completely clumpy (lightest gray; red in color version) with nearly evacuated voids to completely 
smooth (darkest grey; dark blue in color version). Each colored line in the upper plot represents 
the average of 200 viewing angles, and the gray scale shows the range spanned by those individual 
sightlines (as the average sightline is not observable from earth). In the upper plot the gray scale 
is the range spanned by the model with 1% smooth-to-clumpy ratio (very evacuated voids), and in 
the lower plot the gray scale is for the 90% smooth model. Although there are still large variations 
with viewing angle, the smoothness of the cocoon is a large enough effect to be observationally 
probed. 
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Fig. 14. — The left panel of this figure is similar to the lower right panel of Fig. 10, showing 
244 models with different distributions of the same total mass of clumpy clouds, but here different 
ratios of smooth to clumpy material are shown, 0.03 (blue circles), 0.1 (squares), 0.3 (triangles), 
and 0.9 (red circles). The right panel collapses the horizontal axis, showing histograms of the 
observable quantity, [5.8]- [24] color, labeled according to smooth-to-clumpy ratio. Clearly, the 
color can diagnose the degree of dumpiness of the medium (if the luminosity of the central source 
is known and the total mass can be estimated e.g. from FIR observations). 
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color box is plotted along with the colors of field objects and UCH II regions. The latter 
are from the survey of Kurtz et al. (1994) that includes the WC89 sample as a subset. Our 
models have IRAS colors that fit the WC89 criteria, but smooth models do also. The colors 
are shown along different sightlines for a model with smooth-to-clumpy ratio of 0.03, 0.1, and 
0.9. In this part of parameter space (hot protostars with a circumstellar dust cocoon about 
a parsec in size with mean density a few times 10 4 cm -3 ), making the cocoon smoother or 
increasing its density while maintaining the same size both result in bluer [12]- [25] colors, 
but not much change in [60]- [12]. Decreasing the total mass without changing the average 
density (a smaller envelope) results in bluer [60]- [12] without a large change in [25]- [12]. 

Spectroscopy can almost always constrain models more tightly than photometry, so we 
next compare against ground-based MIR spectroscopy of UCH II regions by F98. We begin 
with G5. 89-0. 39, an UCH II which was given special attention by those authors and which 
had been modeled previously. Figure 16 shows data from F98, IRAS, and MSX along with 
a smooth model and our canonical clumpy model. We compared the data with our grid 
of smooth models and confirm that the parameters similar to those given by F98 (they 
used rj n =0.032pc, r out =2.6pc, a=0.0, M=70000M Q ) produce the best fit, but that fit is 
still not perfect, in particular shortward of 8/xm and around 20/im. We then ran a large 
number (~600) of clumpy models with different clump distributions, inner radii, average 
radial density gradients, and smooth-to-clumpy ratios (but fixed outer radius), and found 
the best-fitting sightline by least-squared difference between the data and model. We fit 
only the F98 data between 3 and ISfim, but the best-fitting models also fit the 
longer wavelength data, with no further adjustment. In particular, the slope of 
the SED around 20/im and the behavior shortward of 8/xm resemble the data much more 
closely than the smooth model. The best-fitting smooth and clumpy models' parameters are 
listed in Table 2. The figure also shows the SEDs for other sightlines of the same clumpy 
model as grayscale, showing that if this is indeed an accurate model of G5.89, that it would 
look rather different from a different vantage point. Figure 17 shows three other objects 
from F98 with different types of SEDs (NIR intensity and silicate emission or absorption), 
all of which are fit better by clumpy models than smooth ones (see x 2 values and models 
parameters in Table 2). In some cases, the model flux falls short of the observations at longer 
wavelengths, in particular at 100/xm. At the size of the 100/xm IRAS beam, the Galaxy is 
very confused, and nearby sources almost certainly contribute some heating of the outer dust 
cocoon, increasing the far-IR flux. 

Figures 18 and 19 show contours of goodness-of-fit for the smooth and clumpy models, 
respectively. Smooth models have the same difficulties fitting real data that they did in 
fitting SEDs from 3-D models above: The better-fitting models have flat density profiles or 
even ones that increase as a function of radius. The outer radii are uncomfortably large, and 
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Fig. 15. — IRAS colors of models and real objects in the Galaxy. This color-color diagram is similar 
to that in WC89, and shows a random sample of field objects (plus signs) along with UHH II regions 
from Kurtz et al. (1994) and WC89 (diamonds). Colors of a model viewed from different angles are 
shown for smooth-to-clumpy ratio of 0.03 (triangles), 0.1 (circles), and 0.9 (squares). Our models 
have the WC89 colors of UCH II models, but so do smooth models - the IRAS bandpasses are not 
extremely sensitive to dust geometry, although there are clear trends in this color space with the 
smooth-to-clumpy ratio. 



Table 2. Parameters of 1-D and 3-D models fitting data 



object 




model 


density 
index a 


[PC] 


[PC] 


smooth-to 
clumpy 


mass 

M 


Average Ay 
to star 


x 2 


G5. 89-0. 39 


smooth 


unconstrained 


+1.5 


0.00039 


0.39 


1 


1700 


109 


2.7 




smooth 


a < 0, r out < 7 


0.0 


0.0059 


5.9 


1 


2.3xl0 5 


109 


6.1 




clumpy 




0.0 


0.0001 


2.5 


0.3 


50000 


194 


2.5 


G29.96-0.02 


smooth 


unconstrained 


+0.5 


0.0039 


1.1 


1 


3500 


50 


5.8 




smooth 


ol < 0, r out < 7 


0.0 


0.0039 


3.9 


1 


51000 


50 


7.1 




clumpy 




0.0 


0.0001 


2.5 


0.03 


50000 


178 


1.1 


G34.26+0.f5 


smooth 


unconstrained 


-0.5 


0.0039 


3.8 


1 


1.3xl0 5 


150 


11. 




clumpy 




-1.0 


0.003 


2.5 


0.1 


50000 


259 


2.5 


G75.78+0.34 


smooth 


unconstrained 


+0.5 


0.0019 


19. 


1 


2.1xl0 7 


100 


4.7 




smooth 


a < 0, r out < 7 


0.0 


0.0039 


3.9 


1 


2.0xl0 5 


200 


16. 




clumpy 




-0.5 


0.0001 


2.5 


0.3 


50000 


238 


1.7 
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even the best-fitting models have x 2 ?t 7. 

Clumpy models with a smooth-to-clumpy ratio of <50% are preferred, with a loose 
preference for the range 0.1-0.3. Models with highly evacuated voids (3% smooth) have 
more difficulty reproducing the 5/im flux and short-wavelength shape of the SED seen in 
some objects. On the other hand, models too smooth begin to show the shortcomings of the 
completely smooth models and also fit less well. When considering the radiation intensity 
in the interior of an externally lit cloud, Bethell et al (2004) found that a smooth-to-clumpy 
ratio of 0.333 in a hierarchically clumped model produced a similar internal intensity to an 
externally lit turbulent hydrodynamic simulation - corroborating evidence of a low smooth 
fraction. 

Finally, we note that models with fairly flat angle-averaged radial density profiles often 
fit better than those with a strong radial gradient, but that the clumpy models with negative 
gradients (decreasing density with radius) fit better than smooth models with negative gra- 
dients. Both findings are important: some UCH II regions such as G5.89 are well-described 
by a star in a clumpy molecular cloud, without any particularly concentrated envelope. This 
is consistent with the notion that UCH lis are relatively evolved protostars, perhaps with 
accretion already halted (Churchwell 2002). 1 However, radial density profiles have been 
invoked to explain the submm spatial intensity distribution of high-mass protostellar objects 
(e.g. Hatchell et al. 2000). A full analysis of submm radial intensity profiles is beyond the 
scope of this paper, but preliminary investigation finds that clumpy models with a between 
-0.5 and -1.0 fit the "non-peaked" SCUBA sources of (Hatchell et al. 2000) as well as 1-D 
models, and clumpy models with a ~-1.5 fit the "peaked" sources very well, without needing 
to invoke extra emission in a "core" as was required by those authors with 1-D models. 

In Figure 20, we show something even more remarkable: the single clumpy model that 
best fits G5.89 can in fact do a reasonable job fitting F98's entire collection of UCH II regions! 
We normalized each region's data to the 100/xm point, mostly to correct for distance, but this 
also corrects for the difference in spectral type between the regions. For a deeply embedded 
object, the shape of the output (observable) SED shows very little variation whether the 
cocoon is heated by an 07 star or a B2 star, because the NIR/MIR are dominated by 
reprocessed emission in both cases (see discussion in Whitney et al. 2004b). The total 
luminosity of the object naturally scales with the bolometric luminosity of the illuminating 
star. It is very interesting that all of these young protostars are consistent with the same 



1 We note that G5. 89-0. 39 is a multiple source (Feldt et al. 2003; Sollins et al. 2004), which is consistent 
with our interpretation as an illuminating source embedded in a molecular cloud without any particular 
radial density gradient. We are currently investigating models with multiple sources, and although beyond 
the scope of this paper, the SEDs of small clusters show many similar effects to those discussed here. 
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Fig. 16. — Comparison of clumpy models to UCH II region G5. 89-0. 39. Observations are spectra 
from F98 (plus signs) and photometry from IRAS (circles), MSX (squares), and Chini et al. (1986) 
(diamond) . The best-fit sightline of a clumpy model is a solid line and the range of SEDs for different 
sightlines of that same clumpy model (gray scale, darker reflecting a higher density of similar SEDs). 
Also shown is the best-fitting smooth model (dashed line). The best clumpy sightline fits the data 
very well; note that we only tried to fit 1< A <15/im, and that the agreement at longer wavelengths 
occurred automatically. In particular the clumpy model fits better at 1.3mm datapoint. The clumpy 
model is constructed with a varying aperture size appropriate to the instruments being shown. This 
mostly affects the longest and shortest wavelengths. 
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Fig. 17.— Comparison of clumpy models to UCH II regions G29.96-0.02, G34.26+0.15, and 
G75. 78+0.34, annotations as Fig. 16. 
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Fig. 18. — Contours of goodness-of-fit of 1-D smooth models to four UCH II regions. The models 
favor large outer radii and often often large inner radii as well - see text for discussion. 



Fig. 19. — Contours of goodness-of-fit of 3-D clumpy models to four UCH II regions. The models 
favor relatively flat radial profiles and low smooth-to-clumpy ratios, and best fits tend to have lower 
X 2 than the best fitting 1-D models - see text for discussion. 
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model seen from different directions. Observational tests must be done using ensembles of 
objects, and some care should be taken not to overinterpret differences between individual 
objects as differences in evolutionary state or cocoon mass. 
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Fig. 20. — The same clumpy model that fits G5.89 (Fig. 16, top panel) can in fact fit F98's entire 
sample of UCH II regions reasonably well! Data for several regions are here plotted as connected 
plus signs on top of the gray-scale range of sightlines for a single clumpy model. 



6. Conclusions 

We use 3-D radiative transfer models to examine the importance of an inhomogeneous 
circumstellar envelope on the observed properties of young embedded massive stars. The 
spectral energy distributions of 3-D models vary significantly with viewing angle, and differ 
significantly from 1-D models. Some of our most important findings are as follows: 

• The SEDs of embedded stars and young stellar objects are significantly affected by 
clumping of the dust at wavelengths <100/zm. Longer- wavelength emission is primarily 
determined by the total mass of thermally emitting dust rather than its distribution, 
and therefore less affected by clumping. 

• Silicate features (the strongest and most noticeable at 10/xm) can be seen in emission 
or absorption in same object, depending on viewing angle. 
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• Attempts to fit 3-D clumpy objects with 1-D models can mis-estimate the size of the 
region and total mass by more than an order of magnitude. The derived values of the 
region size, mass, central object luminosity, and optical depth are very often incorrect 
by a factor of 2-3. 

• The typical or average SED (considering all viewing angles) is sensitive to the overall 
smoothness of the dust distribution, or the smooth-to-clumpy ratio. The UCH lis 
observed by F98 are best fit by models that are <50% smooth. 

• Our 3-D models do a better job at fitting the entire SED of UCH lis than previous 1-D 
models, and can easily explain the variation in UCH II SEDs as differences in viewing 
angle of a clumpy medium. 

• Using 3-D clumpy models tends to allow the data to be fit well with more peaked 
density profiles than the best 1-D fits, which may agree better with SCUBA imaging 
data. Some UCH II are best fit with flat density profiles, whether one uses 1-D or 
3-D models - these may be more evolved objects best understood as just a hot star 
embedded in a non-accreting cloud. 

It is becoming clear that models with more than one dimension are necessary to under- 
stand the observed properties of interstellar and circumstellar dust and gas. We have shown 
here that dumpiness can severely affect the emergent spectral energy distribution of embed- 
ded massive stars and UCH II regions. Previous work by us and others has shown that at 
least two-dimensional models are critical to understanding accreting protostars; dumpiness 
in such YSOs may have further effects. In the general ISM, dumpiness affects the radial 
profiles and line ratios observed in H II regions Wood et al. (2005). Clumping has been 
invoked to explain the observed silicate feature in AGN torii Nenkova, Ivezic, & Elitzur 
(2005). Models such as our code are invaluable to determine which physical parameters are 
observationally accessible in complex and multi-phase systems. 
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